home *** CD-ROM | disk | FTP | other *** search
/ Languguage OS 2 / Languguage OS II Version 10-94 (Knowledge Media)(1994).ISO / language / parallax / more_exa.tar / more / nbody.p < prev    next >
Text File  |  1991-11-26  |  4KB  |  156 lines
  1. SYSTEM N_BODY_PROBLEM;
  2.  
  3. CONST n        = 3;
  4.       G        = 2.4024E-4; (* km**3 / (Mt * h**2) *)
  5.       hours    = 24;
  6.       zunit    = 100;
  7.       iter     = hours * zunit;
  8.       zstep    = 1.0 / FLOAT(zunit); (* h *)
  9.       zq_2     = zstep * zstep / 2.0; (* h**2 *)
  10.  
  11. TYPE data = ARRAY [1..n] OF REAL;
  12.  
  13. CONST minit  = data(1.994E21, 5.97E15, 7.35E13);
  14.              (* mass : Mt *)
  15.       xinit  = data(0.0, 1.495E8, 1.491E8);
  16.              (* X-coordinates : km *)
  17.       yinit  = data(0.0, 0.0, 0.0);
  18.              (* Y-coordinates : km *)
  19.       zinit  = data(0.0, 0.0, 0.0);
  20.              (* Z-coordinates : km *)
  21.       vxinit = data(0.0, 0.0, 0.0);
  22.              (* X-speed : km/h *)
  23.       vyinit = data(0.0, 1.072E5, 1.038E5);
  24.              (* Y-speed : km/h *)
  25.       vzinit = data(0.0, 0.0, 2.97E2);
  26.              (* Z-speed : km/h *)
  27.  
  28. SCALAR dvar : data;
  29.        time : INTEGER; (* Time in Hours  * zunit *)
  30.  
  31. CONFIGURATION BODY[1..n];
  32.  
  33. VECTOR mass       : REAL;     (* mass : Mt *)
  34.        x,  y,  z  : REAL;     (* coordinates : km *)
  35.        vx, vy, vz : REAL;     (* speed : km/h *)
  36.        ax, ay, az : REAL;     (* accelerations : km/h**2 *)
  37.  
  38. CONFIGURATION interaction[1..n],[1..n];
  39. CONNECTION BODY1: interaction [i,j] <-> BODY[i].ww1;
  40.            BODY2: interaction [i,j] <-> BODY[j].ww2;
  41.  
  42. VECTOR m1, m2     : REAL;     (* mass : Mt *)
  43.        gm         : REAL;     (* G * m2 : km**3 / h**2 *)
  44.        
  45. SCALAR i : INTEGER;
  46.  
  47. PROCEDURE init;
  48. BEGIN
  49.   time := 0;
  50.   PARALLEL BODY
  51.     dvar := minit;
  52.     LOAD(mass, dvar);
  53.     dvar := xinit;
  54.     LOAD(x    , dvar);
  55.     dvar := yinit;
  56.     LOAD(y    , dvar);
  57.     dvar := zinit;
  58.     LOAD(z    , dvar);
  59.     dvar := vxinit;
  60.     LOAD(vx   , dvar);
  61.     dvar := vyinit;
  62.     LOAD(vy   , dvar);
  63.     dvar := vzinit;
  64.     LOAD(vz   , dvar);
  65.     SEND BODY.ww1 (mass) TO interaction.BODY1 (m1);    
  66.     SEND BODY.ww2 (mass) TO interaction.BODY2 (m2);    
  67.   ENDPARALLEL;
  68.  
  69.   PARALLEL interaction
  70.     gm := G * m2;
  71.   ENDPARALLEL
  72. END init;
  73.  
  74. PROCEDURE step;
  75.  
  76. VECTOR x1, y1, z1 : REAL;     (* coordinates BODY1 : km *)
  77.        x2, y2, z2 : REAL;     (* coordinates BODY2 : km *)
  78.        dx, dy, dz : REAL;     (* distance-components : km *)
  79.        dq         : REAL;     (* square of distance : km**2 *)
  80.        k          : REAL;     (* common constants : 1/h**2 *)
  81.  
  82. BEGIN
  83.   INC (time);
  84.   PARALLEL interaction[*],[DIM1 <> DIM2]
  85.     RECEIVE interaction.BODY1 (x1) FROM BODY.ww1 (x);
  86.     RECEIVE interaction.BODY2 (x2) FROM BODY.ww2 (x);
  87.     RECEIVE interaction.BODY1 (y1) FROM BODY.ww1 (y);
  88.     RECEIVE interaction.BODY2 (y2) FROM BODY.ww2 (y);
  89.     RECEIVE interaction.BODY1 (z1) FROM BODY.ww1 (z);
  90.     RECEIVE interaction.BODY2 (z2) FROM BODY.ww2 (z);
  91.     dx := x2 - x1;
  92.     dy := y2 - y1;
  93.     dz := z2 - z1;
  94.     dq := dx*dx + dy*dy + dz*dz;
  95.     k := gm / (dq * Sqrt(dq));
  96.     SEND interaction.BODY1 (k * dx) TO BODY.ww1 (ax) REDUCE.SUM;
  97.     SEND interaction.BODY1 (k * dy) TO BODY.ww1 (ay) REDUCE.SUM;
  98.     SEND interaction.BODY1 (k * dz) TO BODY.ww1 (az) REDUCE.SUM;
  99.   ENDPARALLEL;
  100.  
  101.   PARALLEL BODY
  102.     x := x + vx * zstep + ax * zq_2;
  103.     y := y + vy * zstep + ay * zq_2;
  104.     z := z + vz * zstep + az * zq_2;
  105.     vx := vx + ax * zstep;
  106.     vy := vy + ay * zstep;
  107.     vz := vz + az * zstep;
  108.   ENDPARALLEL
  109. END step;
  110.  
  111. PROCEDURE output;
  112.  
  113. TYPE name = ARRAY [1..8] OF CHAR;
  114.  
  115. SCALAR dfeld : data;
  116.        i : INTEGER;
  117.  
  118.   PROCEDURE outline(SCALAR s : name);
  119.   BEGIN
  120.     WriteString(s);
  121.     FOR i := 1 TO n DO
  122.       WriteString("   ");
  123.       WriteReal(dfeld[i],15);
  124.     END; (* FOR *)
  125.     WriteLn;
  126.   END outline;
  127.   
  128. BEGIN
  129.   WriteInt(time DIV zunit,5);
  130.   WriteString("    ");
  131.   FOR i := 1 TO n DO
  132.     WriteString("   BODY"); WriteInt(i,3); WriteString("     ");
  133.   END; (* FOR *)
  134.   WriteLn;
  135.   PARALLEL BODY
  136.     STORE(mass , dfeld); outline("m  :  ");
  137.     STORE(x    , dfeld); outline("x  :  ");
  138.     STORE(y    , dfeld); outline("y  :  ");
  139.     STORE(z    , dfeld); outline("z  :  ");
  140.     STORE(vx   , dfeld); outline("vx :  ");
  141.     STORE(vy   , dfeld); outline("vy :  ");
  142.     STORE(vz   , dfeld); outline("vz :  ");
  143.   ENDPARALLEL
  144. END output;
  145.  
  146. BEGIN
  147.   init;
  148.   output;
  149.   LOOP
  150.     FOR i := 1 TO iter DO
  151.       step;
  152.     END; (* FOR *) 
  153.     output;
  154.   END; (* LOOP *)
  155. END N_BODY_PROBLEM.
  156.